home *** CD-ROM | disk | FTP | other *** search
/ Aminet 33 / Aminet 33 - October 1999.iso / Aminet / dev / c / GAPLib.lha / GAPLib / diagnostic / GaussRandTest.c < prev    next >
Encoding:
C/C++ Source or Header  |  1999-05-23  |  1.1 KB  |  71 lines
  1. #include <stdio.h>
  2. #include <stdlib.h>
  3. #include <math.h>
  4. #include <time.h>
  5. #include <GAP.h>
  6.  
  7. #define    LOOPS 1048576
  8.  
  9. int main(void)
  10. {
  11. int ret=0;
  12. int i,n,t;
  13. int    *sumar;
  14. double d,q,err;
  15.  
  16. InitRand(time(0));
  17.  
  18. sumar = malloc(sizeof(double)*129);
  19.  
  20.  
  21. if(sumar!=0) {
  22.  
  23.     t=0;
  24.     for(i=0;i!=129;i++) {
  25.         sumar[i] = 0.0;
  26.     }
  27.  
  28.     for(i=0;i!=LOOPS;i++) {
  29.         d = GaussRand(0.0,0.8);
  30.         n = (int)(32.0*d)+64;
  31.         if(n>=0 && n<128) {
  32.             sumar[n] = sumar[n] + 1.0;
  33.         }
  34.         if((i&4095)==4095) {
  35.             printf("%c\b","-\\|/-\\|/"[t++]);
  36.             fflush(stdout);
  37.             t&=7;
  38.         }
  39.     }
  40.  
  41.     /*
  42.     *    It seems that the average value of the gaussian noise appears
  43.     * exactly twice as often as it should.
  44.     */
  45.  
  46.     sumar[64] = sumar[64]/2.0;    /* This feels like a bug. */
  47.  
  48.     q=1.0/(sumar[64]);
  49.     err = 0.0;
  50.  
  51.     for(i=0;i!=129;i++) {
  52.         d = exp(-pow(i/32.0-2,2.0));
  53.         err = err + pow((d-sumar[i]*q),2.0);
  54.     }
  55.  
  56.     printf("(%f)\n",err); 
  57.     if(err<1.41421356) {
  58.         printf("GaussRandTest: No worse than usual.\n");
  59.     } else {
  60.         printf("GaussRandTest: This sucks.\n");
  61.         ret = 5;
  62.     }
  63.  
  64.     free(sumar);
  65. } else {
  66.     fprintf(stderr,"***Error: No free store. (malloc() failed)\n");
  67. }
  68.  
  69. return(ret);
  70. }
  71.